#include "Splines.h"
#include <fstream>
#include <string>

int main(int argc, char *argv[])
{
    Complete_cubic_spline<double> CX;
    Complete_cubic_spline<double> CY;
    int N;
    std::vector<double> x;
    std::vector<double> y;
    std::vector<double> t;
    double xdfa;
    double xdfb;
    double ydfa;
    double ydfb;
    std::vector<double> pt(500);
    std::string tmp;

    //curve1
    std::ifstream data1("data/curve1.txt");
    data1 >> N;
    std::getline(data1, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data1 >> x[i];
	data1 >> y[i];
	std::getline(data1, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px1(500);
    std::vector<double> py1(500);
    for(int i = 0; i < 500; i++)
    {
	px1[i] = CX.Spline_generator(pt[i]);
	py1[i] = CY.Spline_generator(pt[i]);
    }

    //curve2
    std::ifstream data2("data/curve2.txt");
    data2 >> N;
    std::getline(data2, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data2 >> x[i];
	data2 >> y[i];
	std::getline(data2, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px2(500);
    std::vector<double> py2(500);
    for(int i = 0; i < 500; i++)
    {
	px2[i] = CX.Spline_generator(pt[i]);
	py2[i] = CY.Spline_generator(pt[i]);
    }

    //curve3
    std::ifstream data3("data/curve3.txt");
    data3 >> N;
    std::getline(data3, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data3 >> x[i];
	data3 >> y[i];
	std::getline(data3, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px3(500);
    std::vector<double> py3(500);
    for(int i = 0; i < 500; i++)
    {
	px3[i] = CX.Spline_generator(pt[i]);
	py3[i] = CY.Spline_generator(pt[i]);
    }

    //curve4
    std::ifstream data4("data/curve4.txt");
    data4 >> N;
    std::getline(data4, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data4 >> x[i];
	data4 >> y[i];
	std::getline(data4, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px4(500);
    std::vector<double> py4(500);
    for(int i = 0; i < 500; i++)
    {
	px4[i] = CX.Spline_generator(pt[i]);
	py4[i] = CY.Spline_generator(pt[i]);
    }

    //curve5
    std::ifstream data5("data/curve5.txt");
    data5 >> N;
    std::getline(data5, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data5 >> x[i];
	data5 >> y[i];
	std::getline(data5, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px5(500);
    std::vector<double> py5(500);
    for(int i = 0; i < 500; i++)
    {
	px5[i] = CX.Spline_generator(pt[i]);
	py5[i] = CY.Spline_generator(pt[i]);
    }

    //curve6
    std::ifstream data6("data/curve6.txt");
    data6 >> N;
    std::getline(data6, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data6 >> x[i];
	data6 >> y[i];
	std::getline(data6, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px6(500);
    std::vector<double> py6(500);
    for(int i = 0; i < 500; i++)
    {
	px6[i] = CX.Spline_generator(pt[i]);
	py6[i] = CY.Spline_generator(pt[i]);
    }

    //curve7
    std::ifstream data7("data/curve7.txt");
    data7 >> N;
    std::getline(data7, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data7 >> x[i];
	data7 >> y[i];
	std::getline(data7, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px7(500);
    std::vector<double> py7(500);
    for(int i = 0; i < 500; i++)
    {
	px7[i] = CX.Spline_generator(pt[i]);
	py7[i] = CY.Spline_generator(pt[i]);
    }

    //curve8
    std::ifstream data8("data/curve8.txt");
    data8 >> N;
    std::getline(data8, tmp);
    CX.set_N(N+1);
    x.resize(N+1);
    CY.set_N(N+1);
    y.resize(N+1);
    t.resize(N+1);
    for(int i = 0; i < N; i++)
    {
	data8 >> x[i];
	data8 >> y[i];
	std::getline(data8, tmp);
    }
    x[N] = x[0];
    y[N] = y[0];
    for(int i = 0; i < N+1; i++)
    {
	if(i == 0)
	    t[i] = 0;
	else
	    t[i] = t[i-1] + sqrt((x[i]-x[i-1])*(x[i]-x[i-1]) + (y[i]-y[i-1])*(y[i]-y[i-1]));
    }
    xdfa = (x[1]-x[0])/(t[1]-t[0]);
    xdfb = (x[N+1]-x[N])/(t[N+1]-t[N]);
    ydfa = (y[1]-y[0])/(t[1]-t[0]);
    ydfb = (y[N+1]-y[N])/(t[N+1]-t[N]);
    CX.set_knots(t);
    CX.set_boundary_condition(xdfa, xdfb);
    CX.set_interpolation_condition(x);
    CX.set_derivatives();
    CX.set_lambda();
    CX.set_mu();
    CY.set_knots(t);
    CY.set_boundary_condition(ydfa, ydfb);
    CY.set_interpolation_condition(y);
    CY.set_derivatives();
    CY.set_lambda();
    CY.set_mu();
    
    for(int i = 0; i < 500; i++)
	pt[i] = t[0] + i*(t[N] - t[0])/499.0;
    std::vector<double> px8(500);
    std::vector<double> py8(500);
    for(int i = 0; i < 500; i++)
    {
	px8[i] = CX.Spline_generator(pt[i]);
	py8[i] = CY.Spline_generator(pt[i]);
    }
    
    std::ofstream output("panda.m");
    
    output << "px1 = [" ;
    for(int i = 0; i < 500; i++)
	output << px1[i] << " ";
    output << "];" << std::endl;
    output << "py1 = [";
    for(int i = 0; i < 500; i++)
	output << py1[i] << " ";
    output << "];" << std::endl;
    output << "plot(px1, py1,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px2 = [" ;
    for(int i = 0; i < 500; i++)
	output << px2[i] << " ";
    output << "];" << std::endl;
    output << "py2 = [";
    for(int i = 0; i < 500; i++)
	output << py2[i] << " ";
    output << "];" << std::endl;
    output << "plot(px2, py2,'k');" << std::endl;
    output << "hold on" << std::endl;
    
    output << "px3 = [" ;
    for(int i = 0; i < 500; i++)
	output << px3[i] << " ";
    output << "];" << std::endl;
    output << "py3 = [";
    for(int i = 0; i < 500; i++)
	output << py3[i] << " ";
    output << "];" << std::endl;
    output << "plot(px3, py3,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px4 = [" ;
    for(int i = 0; i < 500; i++)
	output << px4[i] << " ";
    output << "];" << std::endl;
    output << "py4 = [";
    for(int i = 0; i < 500; i++)
	output << py4[i] << " ";
    output << "];" << std::endl;
    output << "plot(px4, py4,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px5 = [" ;
    for(int i = 0; i < 500; i++)
	output << px5[i] << " ";
    output << "];" << std::endl;
    output << "py5 = [";
    for(int i = 0; i < 500; i++)
	output << py5[i] << " ";
    output << "];" << std::endl;
    output << "plot(px5, py5,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px6 = [" ;
    for(int i = 0; i < 500; i++)
	output << px6[i] << " ";
    output << "];" << std::endl;
    output << "py6 = [";
    for(int i = 0; i < 500; i++)
	output << py6[i] << " ";
    output << "];" << std::endl;
    output << "plot(px6, py6,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px7 = [" ;
    for(int i = 0; i < 500; i++)
	output << px7[i] << " ";
    output << "];" << std::endl;
    output << "py7 = [";
    for(int i = 0; i < 500; i++)
	output << py7[i] << " ";
    output << "];" << std::endl;
    output << "plot(px7, py7,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px8 = [" ;
    for(int i = 0; i < 500; i++)
	output << px8[i] << " ";
    output << "];" << std::endl;
    output << "py8 = [";
    for(int i = 0; i < 500; i++)
	output << py8[i] << " ";
    output << "];" << std::endl;
    output << "plot(px8, py8,'k');" << std::endl;
    output << "hold on" << std::endl;

    output << "px = zeros(1,500);" << std::endl;
    output << "py = zeros(1,500);" << std::endl;
    
    output << "k = 1;" << std::endl;
    output << "for i=0:length(px4)-1" << std::endl;
    output << "px(k) = px4(length(px4)-i);" << std::endl;
    output << "py(k) = py4(length(px4)-i);" << std::endl;
    output << "k = k+1;" << std::endl;
    output << "end" << std::endl;
    output << "fill([px3,px],[py3,py],'k');" << std::endl;
    
    output << "k = 1;" << std::endl;
    output << "for i=0:length(px5)-1" << std::endl;
    output << "px(k) = px5(length(px5)-i);" << std::endl;
    output << "py(k) = py5(length(px5)-i);" << std::endl;
    output << "k = k+1;" << std::endl;
    output << "end" << std::endl;
    output << "fill([px1,px],[py1,py],'k');" << std::endl;

    output << "k = 1;" << std::endl;
    output << "for i=0:length(px6)-1" << std::endl;
    output << "px(k) = px6(length(px6)-i);" << std::endl;
    output << "py(k) = py6(length(px6)-i);" << std::endl;
    output << "k = k+1;" << std::endl;
    output << "end" << std::endl;
    output << "fill([px2,px],[py2,py],'k');" << std::endl;

    output << "k = 1;" << std::endl;
    output << "for i=0:length(px8)-1" << std::endl;
    output << "px(k) = px8(length(px8)-i);" << std::endl;
    output << "py(k) = py8(length(px8)-i);" << std::endl;
    output << "k = k+1;" << std::endl;
    output << "end" << std::endl;
    output << "fill([px7,px],[py7,py],'k');" << std::endl;
    
    output.close();

    return 0;
};
